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Abstract-  This  paper  presents  the  development  of  a  signal 
averaging  algorithm  for  recovering  excitation  responses 
contaminated  by  overwhelming  amount  of  various  types 
of  interference  in  skin  admittance  measurements.  The 
algorithm  is  designed  to  eliminate  Gaussian-distributed 
noise  by  use  of  a  recursive  approach.  The  process  of 
recovering  low  magnitude  voltage  responses  from  highly 
noise-contaminated  waveforms  is  a  CPU-intensive  task.  In 
real-time  measurements,  iterative  reconstruction 
algorithm  is  ineffident  and  time  consuming  when  slow 
varying  input  waveforms  are  present.  To  increase  the 
quality  of  the  reconstruction  a  considerably  large  number 
of  recursions  is  required.  Increasing  the  number  of 
recursions  is  appropriate  for  batch  processing  of 
measurement  data.  However,  the  algorithm  considers 
measurements  in  real-time,  whereas  required  quality  of 
signal  reconstruction  should  be  kept  independent  from 
the  number  of  recursions . 

I.  INTRODUCTION 

For  medical  measurement  systems  choosing  appropriate  filter 
specifications  is  not  always  a  straightforward  process.  The 
difficulty  generally  arises  from  insufficient  knowledge  of  the 
unknown  signal  power  spectrum  and  the  noise  power 
spectrum.  The  major  difficulty  arises  when  noise  signal  is 

randomly  occurring  with  a  wide  frequency  range  that 
overwhelms  the  frequency  range  of  the  measured  signal.  In 
this  case  the  type  of  filter  implementation  is  a  dominant  part 
of  the  overall  system  design.  Signal  averaging  provides 
excellent  results  for  noise  removal,  when  the  signal  is 

corrupted  by  additive  white  Gaussian  noise.  Recently,  there 

have  been  many  attempts  to  improve  the  denoising 
performance  at  small  sample  sizes  by  using  statistical 
inference  methods  based  on  wavelet  statistical  models  and 
Bayesian  estimation  [11,12].  The  proposed  algorithm  is 
rather  different  and  will  efficiently  increase  the  S/N  ratio  by 
averaging  out  the  unwanted  signal,  even  with  very  low  level 
of  input  signal  amplitudes. 

II.  M  EASUREMENT  SETUP 

Measurements  were  performed  by  simultaneous  recording  of 
dry  skin  parameters.  The  skin  parameters  of  interest,  as 

described  in  [2],  are  D.C  potential,  A.C  conductance, 
capacitance,  and  the  changes  in  these  caused  by  the  reflex. 
A.C  measurements  were  performed  by  applying  low 
amplitude  sinusoidal  in  the  frequency  range  of  0.1-1000  Hz. 
The  measurement  circuit  is  based  on  three-electrode  system 
[7]  with  constant  current,  which  were  designed  to  record 
simultaneous  individual  measurement  of  electrodes.  The 


measurement  responses  at  different  frequencies  were  A/D 
converted  and  recorded  for  denoising  process  by  use  of  the 
averaging  algorithm. 

HI.  METHOD 

Signal  averaging  is  often  done  by  a  dedicated  computer  after 
A/D  conversion.  However,  averaging  requires  that  a  large 
number  of  bits  per  unit  time  be  processed.  This,  in  turn, 
requires  a  fast  A/D  conversion  and  CPU -intensive  digital  data 
processing.  The  measurement  system  was  simulated  with 
several  different  system  configurations,  and  several  optimum 
results  were  achieved.  By  repetitive  additions  of  waveforms, 
random  noise  samples  tend  to  average  to  zero  while  the 
amplitude  of  the  desired  signal  increased  with  decreasing 
system  bandwidth  and  performance.  It  is  assumed  that  the 
random  noise  signal  must  not  be  correlated  with  the  desired 
signal.  Overlapping  frequency  components  are  not  cut  off,  on 
the  contrary,  they  are  manipulated  equally,  i.e.,  increased  in 
amplitude  along  with  the  input  signal. 

In  order  to  eliminate  the  random  noise  samples  efficiently,  it 
was  necessary  to  increase  the  number  of  additions.  This 
reduces  the  system  bandwidth  unduly  and  causes  unstable 
conditions.  To  overcome  such  situations,  additional 
operations  have  been  applied  in  parallel  with  the  averaging 
process.  Fig.  1  illustrates  an  exUeme  situation,  where  a  256- 
samples  waveform  with  a  S/N  ratio  of  1/100  is  given  as  the 
input. 

Noisy  Input,  S/N=1/100,Samples=256 


Samples 


Fig.  1  Noise  corrupted  waveform  with  256  samples  and  S/N  ratio  of  1/100. 
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The  output  waveform  after  averaging  process  is  shown  in  Fig. 
2,  which  shows  that  the  averaging  has  its  superficial  result,  in 
which  extreme  spikes  of  random  noise  have  vanished 
satisfactorily  only  after  50.000  recursions. 


Recursions=50000,  Samples=256 


Fig.  2  Averaged  waveform  of  the  noise  connpted  waveform  with  256 
samples.  The  connpted  waveform  is  nearly  reconstructed  after  50.000 
iterations. 

IV.  APPROACH 

Let  T  denote  the  time  for  each  waveform,  and  if  M 
waveforms  of  the  noise-contaminated  signal,  each  of  duration 
T  will  be  averaged,  N  samples  from  each  waveform  must  be 
taken  and  stored,  giving  M  XN  samples  in  total,  then  a 
synchronized  addition  of  the  samples  can  be  performed  to 
produce  the  average  value  of  the  M  waveforms.  The  addition 
of  the  samples  can  be  exactly  sychronised  if  the  t1  sample  of 
the  first  waveform  ( Hi )  is  taken  at  time  1 0 ,  the  1th  sample  of 

the  second  waveform  (w-,)  taken  at  time  t0  +  T ,  and  the  i* 

sample  of  the  bf1  waveform  ( Wn  )  taken  at  time  t0  +  NT  . 
As  shown  in  the  following  Mathematica  code  slice,  each 
waveform  (vt;,  )  is  an  obserx'ed  sample  denoting  a  waveform 
of  a  periodic  signal  containing  N  discrete  random  noise  and 
measured  signal  samples. 

NoiseSamp/:  NoiseSamp[magnJ:=  magn*Random[Real,{-1 ,1}] 
InpSamp/:  lnpSamp[sampno_lnteger,n 

magn_,samps_lnteger,freq_:1  ,Phase_:0]  := 

N[Sin[2  Pi(sampno*freq+Phase)/samps],7]+NoiseSamp[nmagn]] 


Where,  s{  f)  denotes  the  signal  sample,  n  (i)  the  noise  sample 

of  a  random  magnitude,  and  Vt>  denotes  sum  of  these 
quantities. 

Corresponding  samples  from  each  waveform  are  then 
averaged  to  give  the  average  value  wk  at  the  k*  sample 
position  of  M  waveforms 

1  M 

Wk  =  17  12 l^n<N  (1) 

At  r=l 

where  W(l,  r)  represents  the  n*  sample  of  the  rth  waveform.  In 
the  simplest  manner,  the  average  value  is  calculated  by 
adding  the  waveforms  , . . . ,  WM  in  a  sample-by-sample 
manner,  and  then  dividing  the  sum  by  the  number  of 
waveforms: 

l  M 

y^M)  =JfY<W(n,r)  n  =l,...,N 

1  N  M 

^(NM  )  -  T7  52  12  W(n,r)  (2) 

1V1  „=1  r=l 

where 

•  W(„  r)  =  nth  sample  of  the  r*  waveform, 

•  Average  value  of  M  waveforms  containing 
N  samples  each. 

The  algorithm  was  first  implemented  and  tested  by  using  the 
following  Mathematica  code  slice: 

Avr /:  Avr[samples_lnteger,recursions_lnteger, 
NoiseLvl_,Frq_:1,Phase_:0]  := 

Block[{n  =  0,r  =  0,  myi,  iterator  =  0,i,  S, 

FilteredWave,  Signalin},  Array[W, samples]; 
lnPut=Table[lnpSamp[i,  NoiseLvl,  samples,  Frq,  Phase], 

]i,0, samples  -1  }]; 

Signalin  =  Table[lnpSamp[i,0, samples, Frq, Phase],  {i,0, samples  - 1}]; 
For  [r=1 ,  r  <=  recursions,  r++,  For  [n=0,  n  <  samples,  n++, 

W[n]=W[n]+(lnpSamp[n, NoiseLvl, samples, Frq, Phase]  -  W[n])/r]; 
iterator++; 

lf[Mod[iterator,  1 0]  ==  0, 

S=Table[W[myi],  ]myi,0,  samples  - 1}]; 

]]; 

FilteredWave=Table[W[myi],  ]myi,0,  samples  - 1}]; 

Return  [FilteredWave]  ] 

The  M  XN  noisy  samples  are  grouped  into  random 
samples  to  create  a  random  space  of  M  observable 
waveforms  defined  as 


This  function  is  a  symmetrical  implementation  of  a  waveform 
of  which  the  half  period  is  calculated  by, 

N 

^=120(0 +"<o) 

i=i 


Wl  —  x(\,\)1  ■  • 

•  ^(lgV) 

W  -  x 
v,2  a(2,1)’- 

•  '  ’X(2,N) 
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fWx  A 


kWMj 


•*(1,1) ' ' '  -*(i,An 


-*(2,1)  '  "X(2,N) 


yX(M,  1)  ‘  ’  ‘  X(M,N)  j 


(3) 


5.= 


fx  ^ 
Ao,i) 


V(2.1) 


\^X(M  ,1)  j 


-,SN  = 


( x  ^ 

(1JV) 


V(2,AT) 


yX(M,N)  J 


(4) 


Y  M 

Ji  =  — /  x(rl 
M  j*  ('J 

Y  M 

s2=  —  Y  x(r: 

-  m  tr 


will  yield  the  average  value,  and  the  consequent  averages  for 
each  waveform  with  n  samples  each,  i.e.. 


W,  ={  wl,w2,---w„  } 

i  N 

w,.=— Vx„ 


(8) 


or  arranging  the  waveforms  into  N  random  observable 
columns 


Recursive  additions  of  M  waveforms  with  N  discrete  samples 
will  yield 


1  M 

w^—Yw,  , 

M tr 


_  J  M 

n= 1 


Thus,  the  average  value  of  the  random  observable  columns 
can  be  determined  as 


n=l,---,N  (9) 

_  ]  M  _ 


_  1  M 

w  = — VL 

M  W) 

_  i  M 

w  =  -Yx(n  r) 
Mtr 

^  =  w>. 


H - +  X. 


(, n,M ) 


) 


n  =  N 


(10) 


=  J_V 

Mh 


\r.N ) 


(5) 


It  is  obvious  that  the  mean  values  Sx ,  •  •  • ,  SN  construct  the 
sample  mean  of  the  M  waveforms,  Wj ,  •  •  • ,  WM  i.e., 

r={s„s2,-,s„ ) 


w.± 

M 


fYM  w  N 

Lr= 1  VV(r,  1) 

YM  w 

Ljy= 1  ( r,Z ) 


£>< 


(r,iV) 


(6) 


v.  Conclusions 

Sample-by-sample  addition  of  two  discrete-time  signals 
derived  as 

Y(n  ,r)  =  % i(„,r )  +  X2  („,r)  (11) 

where 

•  r)  =  Sum  output  of  the  nlh  sample  at  the  rth 
recursion, 

*  XHnr)=  nth  input  sample  of  the  discrete-time 
signal  >1(77  ]  at  the  Ith  recursion. 

Generally,  for  M  recursions  the  output  sum  is  written  as 


M 

Y(N ,M)  =  ^  (wi(B,7-)  4-  W2(n 


,0)  n  =  !,■■■,  N  (12) 


A  rather  intuitive  alternative  of  determining  the  sample  mean 
W  may  be  acquired  by  taking  average  of  each  waveform  as 
a  whole,  and  measuring  the  mean  value  of  the  averaged 
waveforms.  Thus,  for  the  M  waveforms 


_  1  M 

w  =—Ywr 

m tr 


(V) 


However,  repetitive  additions  of  discrete-time  samples  have 
produced  undesired  output  levels,  which  might  be  difficult  to 
manipulate  by  electronic  circuits.  Increased  number  of 
recursions  can  cause  noticeable  degradation  in  system 
dynamics.  We  have  therefore  introduced  an  appropriate 
expression  in  order  to  keep  the  output  under  control.  When 
the  interference  was  high  and  thousands  of  iterations  were 
required  the  control  mechanism  operated  well.  Therefore,  this 
algorithm  applies  equally  well  to  both  batch  processing  and 
real-time  measurements. 


Proceedings  -  23rd  Annual  Conference  -  IEEE/EMBS  Oct.25-28,  2001,  Istanbul,  TURKEY 


Keeping  this  in  mind,  reasonable  output  levels  were  acquired 
by  modifying  Equation  12.  If,  after  each  addition,  the  output 
were  scaled  down  by  a  certain  factor,  increasing  repetitions 
would  not  cause  any  growth  in  the  output  level.  The  equation 
was  therefore  modified  to  provide  stable  characteristics. 


l(NM) 


=E 


xm-y, 


( n,r-{ ) 


+Y, 


Um-1) 


n  =  !,■■■,  N  (13) 


where  is  substituted  for  the  input  term  w1(„jr)and, 

Y(n  r_tj  is  substituted  for  the  input  term,  W2(n  r)  which  is  the 

sum  of  the  repetitions  prior  to  the  r-lth  recursion.  This 
equation  was  applied  in  the  system  simulation,  which 
provided  efficient  system  dynamics.  Fig.  3  shows  the  result 
of  a  real-time  measurement;  the  noisy  input,  and  the  averaged 
output  after  a  2000  recursions  of  averaging  process. 


Fig.  3  Photo  of  a  real-time  measurement.  Only  noisy  input  and  the  averaged 
output  are  shown. 
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